% Compute Mean-Squared Errors, Bias-Squared, Variance
method_names={'VAR','LP','Pool','Poolfc'};
n_methods=4;

results.irf.(method_names{1})=IRF_VAR_result;
results.irf.(method_names{2})=IRF_LP_result;
results.irf.(method_names{3})=IRF_Pool_result;
results.irf.(method_names{4})=IRF_Pool_fc_result;

the_true_irf = IRF_True';                   % True IRF
the_rms_irf  = sqrt(mean(the_true_irf.^2)); % Root average squared true IRF across horizons
for i_method = 1:n_methods    
    thisMethod = method_names{i_method};    
   results.MSE.(thisMethod) = squeeze(mean((results.irf.(thisMethod) - permute(IRF_True',[1 3 2])).^2, 2));    
   results.BIAS2.(thisMethod) = (squeeze(mean(results.irf.(thisMethod), 2)) - IRF_True').^2;    
   results.VCE.(thisMethod) = squeeze(var(results.irf.(thisMethod), 0, 2));
end

% Criteria results
  RMSE_result= sqrt(extract_struct(results.MSE));
  RMSE_result_plot=squeeze(mean(RMSE_result, 2));
    
  BIAS_result=sqrt(extract_struct(results.BIAS2));
  BIAS_result_plot=squeeze(mean(BIAS_result, 2));
     
  STD_result=sqrt(extract_struct(results.VCE));
  STD_result_plot=squeeze(mean(STD_result, 2));

% Average length
IRF_VAR_seqblock_avglength=mean(IRF_VAR_seqblock_length,2);
IRF_LP_seqblock_avglength=mean(IRF_LP_seqblock_length,2);
IRF_Pool_seqblock_avglength=mean(IRF_Pool_seqblock_length,2);
IRF_Pool_seqblock_avglength_bs1=mean(IRF_Pool_seqblock_length_bs1,2);
IRF_Pool_seqblock_avglength_bs2=mean(IRF_Pool_seqblock_length_bs2,2);

IRF_VAR_bootrep_avglength=mean(IRF_VAR_bootrep_length,2);
IRF_LP_bootrep_avglength=mean(IRF_LP_bootrep_length,2);
IRF_Pool_bootrep_avglength=mean(IRF_Pool_bootrep_length,2);

% Coverage rate
CVG_var_seqblock=zeros(hor+1,n_MC);
CVG_var_bootrep=zeros(hor+1,n_MC);
CVG_lp_seqblock=zeros(hor+1,n_MC);
CVG_lp_bootrep=zeros(hor+1,n_MC);
CVG_pool_seqblock=zeros(hor+1,n_MC);
CVG_pool_seqblock_bs1=zeros(hor+1,n_MC);
CVG_pool_seqblock_bs2=zeros(hor+1,n_MC);
CVG_pool_bootrep=zeros(hor+1,n_MC);

for i_MC = 1:n_MC  
   for j=1:hor+1      
        if  (IRF_True(:,j)>=IRF_VAR_LB_seqblock(j,i_MC)) && (IRF_True(:,j)<=IRF_VAR_UB_seqblock(j,i_MC))
                CVG_var_seqblock(j,i_MC) = 1;
        end

        if  (IRF_True(:,j)>=IRF_VAR_LB_bootrep(j,i_MC)) && (IRF_True(:,j)<=IRF_VAR_UB_bootrep(j,i_MC))
                CVG_var_bootrep(j,i_MC) = 1;
        end

        if  (IRF_True(:,j)>=IRF_LP_LB_seqblock(j,i_MC)) && (IRF_True(:,j)<=IRF_LP_UB_seqblock(j,i_MC))
                CVG_lp_seqblock(j,i_MC) = 1;
        end

        if  (IRF_True(:,j)>=IRF_LP_LB_bootrep(j,i_MC)) && (IRF_True(:,j)<=IRF_LP_UB_bootrep(j,i_MC))
                CVG_lp_bootrep(j,i_MC) = 1;
        end

        if  (IRF_True(:,j)>=IRF_Pool_LB_seqblock(j,i_MC)) && (IRF_True(:,j)<=IRF_Pool_UB_seqblock(j,i_MC))
                CVG_pool_seqblock(j,i_MC) = 1;
        end

        if  (IRF_True(:,j)>=IRF_Pool_LB_seqblock_bs1(j,i_MC)) && (IRF_True(:,j)<=IRF_Pool_UB_seqblock_bs1(j,i_MC))
                CVG_pool_seqblock_bs1(j,i_MC) = 1;
        end

        if  (IRF_True(:,j)>=IRF_Pool_LB_seqblock_bs2(j,i_MC)) && (IRF_True(:,j)<=IRF_Pool_UB_seqblock_bs2(j,i_MC))
                CVG_pool_seqblock_bs2(j,i_MC) = 1;
        end

        if  (IRF_True(:,j)>=IRF_Pool_LB_bootrep(j,i_MC)) && (IRF_True(:,j)<=IRF_Pool_UB_bootrep(j,i_MC))
                CVG_pool_bootrep(j,i_MC) = 1;
        end
                             
   end
end         

% Average coverage
CVG_var_seqblock_avg=sum(CVG_var_seqblock,2)/n_MC;
CVG_var_bootrep_avg=sum(CVG_var_bootrep,2)/n_MC;
CVG_lp_seqblock_avg=sum(CVG_lp_seqblock,2)/n_MC;
CVG_lp_bootrep_avg=sum(CVG_lp_bootrep,2)/n_MC;
CVG_pool_seqblock_avg=sum(CVG_pool_seqblock,2)/n_MC;
CVG_pool_seqblock_avg_bs1=sum(CVG_pool_seqblock_bs1,2)/n_MC;
CVG_pool_seqblock_avg_bs2=sum(CVG_pool_seqblock_bs2,2)/n_MC;
CVG_pool_bootrep_avg=sum(CVG_pool_bootrep,2)/n_MC;


% Average correlation
correlation_varlp_seqblock_save=mean(abs(correlation_varlp_seqblock_sim),4);
correlation_varlp_seqblock_save_bs1=mean(abs(correlation_varlp_seqblock_sim_bs1),4);
correlation_varlp_seqblock_save_bs2=mean(abs(correlation_varlp_seqblock_sim_bs2),4);
correlation_varlp_bootrep_save=mean(abs(correlation_varlp_bootrep_sim),4);